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Abstract 



, Recently, the SPIKE-distance has been proposed as a parameter- free and time-scale independent measure of spike train synchrony. 
"This measure is time-resolved since it relies on instantaneous estimates of spike train dissimilarity. However, its original definition 
led to spuriously high instantaneous values for event-like firing patterns. Here we present a substantial improvement of this measure 
' ^ .which eliminates this shortcoming. The reliability gained allows us to track changes in instantaneous clustering, i.e., time-localized 
^.^patterns of (dis)similarity among multiple spike trains. Additional new features include selective and triggered temporal averaging 
r~] as well as the instantaneous comparison of spike train groups. In a second step, a causal SPIKE-distance is defined such that the 
O ji nstanta.nenus values of dissimilarity rely on past information only so that time-resolved spike train synchrony can be estimated in 
' '"real-time. We demonstrate that these methods are capable of extracting valuable information from field data by monitoring the 
^ synchrony between neuronal spike trains during an epileptic seizure. Finally, the applicability of both the regular and the real-time 
^ SPIKE-distance to continuous data is illustrated on model electroencephalographic (EEG) recordings. 

^: 

^ }Cey words: data analysis; synchronization; spike trains; clustering; SPIKE-distance 
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1. Introduction 

(N 

. . . Measures that estimate the degree of synchrony between 
. ^ y^'^ '^^ more simultaneously recorded spike trains are im- 
portant tools in many different areas of neuroscientific re- 
$^ 'search (Kreuz, 2011). Among many potential applications 
5^ they can be used to evaluate the role of synchronous neu- 
ronal firing in signal propagation (Reyes, 2003), or to esti- 
mate the pairwise correlation within a neuronal population 
(Pillow et al., 2008). 

Some of the most widely used measures of spike train 
dissimilarity depend on a parameter that determines the 
temporal scale in the spike trains to which the measure is 
sensitive. Examples include the Victor-Purpura spike train 
distance (Victor and Purpura, 1996) and the van Rossum 
distance (van Rossum, 2001). Recently, the ISI-distance 
(Kreuz et al., 2007, 2009) and the SPIKE-distance (Kreuz 
et al., 2011) have been proposed as parameter-free and 
time-scale adaptive alternatives. These two measures are 
complementary to each other: While the ISI-distance quan- 
tifies local dissimilarities based on the neurons' firing rate 
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profiles, the SPIKE-distance tracks dyssynchrony mediated 
by differences in spike timing. The latter kind of sensitiv- 
ity can be very relevant since coincident spiking is found in 
many different neuronal circuits, e.g., in the visual cortex 
(Usrey et al., 1999; Priebc and Ferster, 2008) and in the 
retina (Meister and Berry, 1995; Shlens et al., 2008). 

In some situations it is sufficient to evaluate spike train 
synchrony at a rather low temporal resolution, e.g., by 
means of a moving- window analysis where the level of syn- 
chronization within a certain interval is compressed into 
a single value and then compared for successive intervals. 
On the other hand, many applications require a high tem- 
poral resolution, e.g., in order to detect replay of precisely 
timed sequential patterns of neural activity (Ji and Wil- 
son, 2007), to track spike train response variability within 
a neuronal population (Mitchell et al., 2007; Kreuz et al., 
2009), or to understand the role of synchronous firing in the 
neuronal coding of time-dependent stimuli (Miller and Wil- 
son, 2008). In epilepsy research, a high temporal resolution 
could help to gain a deeper understanding of the neuronal 
spiking patterns involved in the different phases of seizure 
generation, propagation, and termination (Truccolo et al., 
2011; Bower et al., 2012). 



The maximum possible temporal resolution is achieved 
when one value of dissimilarity is calculated for each time 
instant and a continuous time profile is obtained. In con- 
trast to the Victor-Purpura and the van Rossum distance, 
the ISI- and the SPIKE-distance can be calculated and vi- 
sualized in such a time-resolved manner. However, as al- 
ready noted in Kreuz et al. (2011), while the original def- 
inition of the SPIKE-distance yields the expected results 
after temporal averaging and also correctly reflects long- 
term trends by means of a moving average, slightly unre- 
liable spiking events lead to spurious high instantaneoiis 
values. In the first part of this study, we remedy this prob- 
lem by considering spike time differences between nearest 
spikes instead of separating differences between preceding 
spikes and differences between following spikes. 

After improving the SPIKE-distance, we extend its ap- 
plicability to situations where the degree of synchrony be- 
tween two or more simultaneously recorded spike trains is 
monitored in real-time. In the field of brain-machine inter- 
faces this may be a promising approach to the rapid on- 
line decoding of neural signals needed to control prosthet- 
ics (Hochberg et al., 2006; Sanchez, 2008). In epileptic pa- 
tients who are refractory to medical treatment, the method 
could be applied to large ensembles of single neurons close 
to or within the epileptic focus and then be integrated into 
a prospective algorithm aimed, e.g., at an early seizure de- 
tection (Jouny et al., 2011). 

Most spike train distances calculate one value of overall 
spike train synchrony for a time interval once this time in- 
terval has passed. There is a lack of measures that can esti- 
mate synchrony in a time-resolved and, even more striking, 
in a causal way. The SPIKE-distance, similar to the ISI- 
distance, is calculated from instantaneous values of spike 
train dissimilarity for which at each time moment not only 
preceding spikes but also following spikes are taken into ac- 
count. This non-causal dependence on future spiking does 
not allow for a real-time calculation. In the second part 
of this study we modify the SPIKE-distance such that the 
instantaneous value of dissimilarity for two or more spike 
trains relies on past information only and can be calculated 
in real-time and in a causal manner. 

We illustrate the new methods on three different appli- 
cations. First we validate the improved definition of the 
SPIKE-distance and its real-time variant on artificially gen- 
erated spike trains. We show that these measures are not 
only able to track the overall synchronization within a 
group of two or more spike trains, but also, due to the re- 
liability added by the revised definition of the distance, to 
track changes in instantaneous clustering, that is, to fol- 
low the evolution of (dis)similarity patterns within multiple 
spike trains. Furthermore, we demonstrate additional fea- 
tures such as selective and triggered temporal averaging as 
well as the comparison of spike train groups. Subsequently, 
we present an application to real spike train data in which 
we analyze single- and multi-unit activity from an epilepsy 
patient recorded in a time interval which includes the oc- 
currence of an epileptic seizure. We can show that the new 



spike train distances are suitable measures to characterize 
the neuronal firing patterns involved in seizure generation, 
propagation, and termination. 

Finally, as for spike trains, there is a lack of measures 
that are capable to monitor time-resolved synchrony in con- 
tinuous data. This issue is addressed in the third and last 
application in which both variants of the SPIKE-distance 
are used to measure the time-resolved dissimilarity in con- 
tinuous data which, in a preprocessing step, are first trans- 
formed into discrete data. We illustrate this approach on 
electroencephalogram (EEG) time series whose level of syn- 
chronization was modified post hoc in a controlled manner 
by means of linear mixing. 

2. Methods 

The different variants of the SPIKE-distance, and the 
ISI-distance (see Appendix A), are defined in terms of a 
function of time, a time profile for each pair of spike trains 
which gives an instantaneous measure of the (dis)similarity 
between the two spike trains. The distances are then defined 
as the temporal average of these dissimilarity profiles, e.g., 
for the bivariate SPIKE-distance, 

Ds=^ r s{t)dt (1) 

^ Jt=o 

where T denotes the overall length of the spike trains, e.g., 
the duration of the recording in an experiment. In the fol- 
lowing this equation is always omitted, and we restrict our- 
selves to showing how to derive the respective dissimilarity 
profiles, e.g., S{t). 

Furthermore, for all distances there exists a straightfor- 
ward extension to the case of more than two spike trains 
(number of spike trains N > 2), the averaged bivariate dis- 
tance. This average over all pairs of neurons commutes with 
the average over time, so it is possible to achieve the same 
kind of time-resolved visualization as in the bivariate case 
by first calculating the instantaneous average, e.g., S'^{t) 
over all pairwise instantaneous values S"^"{t), 

n—l m— n+1 

and only then averaging the resulting dissimilarity profile 
using Eq. 1. All bivariate and averaged bivariate dissim- 
ilarity profiles and thiis all distances arc bounded in the 
interval [0, 1]. The value zero is only obtained for identical 
spike trains. 

2.1. Original definition of the SPIKE-distance 

We first review the original definition of the time- 
resolved profile for the bivariate SPIKE-distance (Kreuz 
et al., 2011). This profile is denoted as So{t) with the 
subscript 'o' standing for original. We then compare it 
against the improved profile, denoted as S{t), by which 
it will henceforth be replaced. The derivation consists of 



three steps: calculating the instantaneous time differences 
between spikes, taking the locally weighted average and 
normalizing the result. 

After labeling the times of the spikes in the spike trains 
n = 1, 2 by i = 1, Af„ (with M„ denoting the num- 
ber of spikes of the n-th spike train), we assign to each time 
instant between and T (Fig. 1) the time of the preceding 
spikes 

(3) 
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the time of the following spikes 



t^^\t)=m.in{tr\tr >t), 

i 

as well as the instantancoiis intcrspike interval 
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The ambiguity regarding the definition of the very first 
and the very last interspike interval as well as the initial 
distance to the preceding spike and the final distance to 
the following spike is resolved by adding auxiliary leading 
spikes at time t = and auxiliary trailing spikes at time 
t = T to each spike train. 

We then denote the instantaneous absolute differences of 
preceding and following spike times as 
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respectively. An instantaneous spike time based measure 
of spike train distance is given by the average of these two 
absolute differences of preceding and following spike times. 
Dividing this average spike time difference by the average of 
the two instantaneous interspike intervals achieves proper 
normalization as well as time-scale invariance (i.e., stretch- 
ing or compressing spike trains does not change the result): 
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In order to be more local in time, a weighted average of 
the two differences Aij(t) with j = P,F is employed such 
that the difference of the spikes that are closer dominate. 
To this aim we denote with 
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and 
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the intervals to the previous and the following spikes for 
each neuron n = 1, 2. 

Inserting the inverse of the mean intervals as weights in 
the locally weighted average 

Ej=P,FAtj(i)- ' 

(Atj(t))j=P,F = 



and making use of 
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Fig. 1. Illustration of the bivariate SPIKE-distance. A. Local quantities 
in relation to time instant t needed to define the original dissimilarity 
profile (Eq. 13) of the SPIKE-distance. B. Additional definitions used 
for the improved dissimilarity profile (Eqs. 17 and 19) and its real-time 
variant (Eq. 21). 



yields the dissimilarity profile 
So{t) = 



Aip(i)(4"Ht))„ + AtF(i)(4"^(t))r 



(Xjgj(f))^ 



(13) 



For multiple spike trains the averaged bivariate variant is 
defined via Eq. 2. 

2.2. Improved definition of the SPIKE-distance 

In its original definition in Kreuz et al. (2011), the 
SPIKE-distance proved to be a reliable indicator of the 
overall level of multi-neuron synchrony for both simulated 
and real data. Appropriate moving averaging also allowed 
to correctly reflect long-term changes in spike train syn- 
chrony. However, as already noted in Kreuz et al. (2011), 
each individual instantaneous value in itself is less reliable 
since spuriously high values are obtained during slightly 
unreliable spiking events. 

This can be seen in the bivariate example of Fig. 2A 
where we use a frequency mismatch to construct two spike 
trains with gradually varying spike matches. In the first half 
the spikes in the second spike train exhibit increasing dis- 
tances to the preceding spikes of the first spike train, while 
in the second half they move closer and closer to its follow- 
ing spikes. Thus we expect a monotonous increase of the 
instantaneous values followed by a monotonous decrease. 
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Fig. 2. Comparison of the original and the improved definition of the 
dissimilarity profile of the SPIKE-distance on constructed spike trains. 
While So{t) shows spurious high values for event-like firing patterns, 
S(t) reflects the level of spike train synchrony faithfully. A. Bivariate 
example: Varying spike matches. B. Multivariate example with 50 spike 
trains: In the first half within the noisy background there are four 
regularly spaced spiking events with increasing jitter. The second half 
consists of ten spiking events with decreasing jitter but now without 
any noisy background. Note that both dissimilarity profiles start at zero 
due to the auxiliary spikes at time t = 0. 



These two trends can indeed be recognized. However they 
are interrupted by short intervals of spurious high values. 

The same kind of spurious high values can be seen in 
Fig. 2B for the averaged bivariate variant Sg{t) in a mul- 
tivariate example. Here in the first half we generated four 
spiking events with increasing jitter within a noisy back- 
ground whereas the second half consists of evenly spaced fir- 
ing events with increasing precision, this time without any 
background noise. In both cases the spurious high values 
appear in small intervals during which some of the spikes 
are still following spikes while others are already preceding 
spikes. In these intervals the small differences between the 
spikes within the event are not taken into account. Instead, 
due to the separation of preceding and following spikes the 
'wrong' spikes are compared to each other. Some of the dif- 
ferences among the preceding spikes as well as among the 
following spikes are very large and this, according to Eq. 
13, leads to the high instantaneous values. 

A straightforward remedy for this mismatch is to com- 
pare the correct spikes to each other, in this case, to com- 
pare each spike to the nearest spike in the other spike train. 
To do so we first rewrite So{t)- For the sake of brevity we 



now omit time dependencies. Using 2(xj"'')„ 
with j = P, F gives 
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and expanding the definition of Atp and Atp (Eqs. 6 and 
7) yields 
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In this formtilation So can be interpreted as the normalized 
sum of four weighted differences, one each for the preced- 
ing spike from the first spike train tp \ the following spike 
from the first spike train tp \ the preceding spike from the 

(2) 

second spike train tp , and, finally, the following spike from 

(2) 

the second spike train tp . However, as we have established 
above, in some instances these four corner spikes are com- 
pared against the 'wrong' spikes. This happens because we 
always restrict the spikes of comparison to the respective 
other preceding and to the respective other following spike. 

A way to resolve this restriction is to allow more flexibil- 
ity and compare each of these four corner spikes to its most 
appropriate spike, i.e., the closest counterpart in the other 
spike train. To this aim we define 
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(16) 



and analogously for tp\ tp\ and tp\ In the improved 
dissimilarity profile these four terms replace the twofold 
contributions of \Atp\ and jAtpl- Furthermore, instead of 
one local weighted average of the two differences between 
previous and following spikes with the mean intervals to the 
previous and the following spikes as weights (Eq. 11), the 
weighting is now carried out for each spike train separately. 
The local weighting for the spike time differences of the first 
spike train reads 
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and analogously 52 (t) is obtained for the second spike train. 

Averaging over the two spike train contributions and nor- 
malizing by the mean interspike interval yields 
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This quantity weights the spike time differences for each 

spike train according to the relative distance of the corner 
spike from the time instant under investigation. This way 
relative distances within each spike train are taken care 
of, while relative distances between spike trains are not. 
In order to get these ratios straight, in a last step the two 
contributions from the two spike trains are locally weighted 
by their instantaneoiis interspike intervals. This leads to 
the improved definition of the SPIKE-distance 

For many time instants the closest spike to the preceding 
spike of one spike train is the preceding spike in the other 
spike train and the same holds for the following spikes, 
so often the same spikes are still compared to each other. 
However, for the time instants that led to the spurious high 
values for So{t) now the two preceding and the two following 
spikes are no longer compared to each other. Instead if a 
following spike in one spike train is the closest spike to a 
preceding spike of the other spike train (or vice versa) these 
will now be compared to each other, which leads to lower 
instantaneous values in the dissimilarity profile. 

This modification resolves the problem of spurious high 
values and, as can be seen in Fig. 2, the desired dis- 
similarity profiles axe obtained. In Fig. 2A the observed 
monotonous increase of the instantaneous values followed 
by a monotonous decrease mirrors exactly the actual 
change in the match between the spike times of the two 
spike trains. For the multivariate example of Fig. 2B and 
the averaged bivariate variant (calculated according to 
Eq. 2) we find in the first half high instantaneous val- 
ues which refiect the background noise. The presence of 
four spiking events with increasing jitter within this noisy 
background is indicated by less and less pronounced drops 
in the dissimilarity profile. In the second half where there 
is no background noise, the evenly spaced firing events 
with increasing precision are correctly reflected by a rather 
monotonous decrease of the instantaneous synchrony level 
which reaches zero at the perfectly synchronous event at 
3900 ms. Such perfect events for which the value zero is 
obtained are marked by vertical dashed lines in the dissim- 
ilarity profile. 

2.3. Real-time SPIKE-distance 

Here we introduce the real-time SPIKE-distance Ds^. 
This is a modification of the SPIKE-distance with the key 
difference that the corresponding time profile Sr {t) can be 
calculated online because it relies on past information only. 
From the perspective of an online measure, the informa- 
tion provided by the following spikes, both their position 
and the length of the interspike interval, is not yet avail- 
able. Like the regular (improved) SPIKE-distance Ds, this 
causal variant is also based on local spike time differences 
but now only two corner spikes are available, and the spikes 



of comparison are restricted to past spikes, e.g., for the pre- 
ceding spike of the first spike train 

a4'^ = min(|4'^ _ t(2)|),t. < t. (20) 

i 

Since there are no following spikes available, there is no 

local weighting, and since there is no interspike interval, the 
normalization is achieved by dividing the average corner 
spike difference by twice the average time interval to the 
preceding spikes (Eq. 9, see also Fig. IB). This yields a 
causal indicator of local spike train dissimilarity: 

S,{t) = . (21) 

In Fig. 3 we show the results of the real-time SPIKE- 
distance for the two cases already used in Fig. 2. As can 
be seen for the example of two spike trains with gradually 
varying spike mismatches (Fig. 3 A) , any spike time differ- 
ence is considered most relevant right at the later of two 
spikes when Sr (t) goes back to a local maximum value. In 
the case where the two preceding spikes are closest to each 
other, it goes back to its maximum value of one. At these 
points the mean time interval to the two preceding spikes is 
exactly half their difference. Any successive period of com- 
mon non-spiking leads to a decrease of the instantaneous 
distance values. This is a desired property since common 
non-spiking is as much a sign of synchrony as common spik- 
ing. The decrease is hyperbolic and its slope depends on 
the preceding spike time difference. 

In addition to the regular trace we here also show a mov- 
ing average which, in line with the real-time calculation, is 
causal. While the regular dissimilarity profile exhibits cer- 
tain fluctuations, this moving average shows that the real- 
time SPIKE-distance in fact reflects the increasing spike 
shifts in the center of the interval and the better match of 
the spike times at the edges. 

During irregular spiking in the multivariate example of 
Fig.3B, the time intervals to the preceding spikes in differ- 
ent spike trains are very variable, and this leads to large 
fluctuating values oi S^{t). Within this irregularity, the pc!r- 
fect spiking event at 400 ms results in an abrupt drop to 
zero, which reflects the delta distribution of the time inter- 
vals to the preceding spikes. The successively more jittered 
events which follow are indicated by a very pronounced 
short-term increase, followed by a decrease. Here the dis- 
tribution of time intervals to preceding spikes starts to de- 
velop a peak at very small intervals and thus becomes bi- 
modal, which causes the increase. Then, after spikes have 
appearc;d in quick succession in all spike trains, the dis- 
tribution becomes very narrow, which is reflected by the 
decrease. In the second half, when there is no background 
noise and the spiking events become less and less noisy, the 
succession of peaks denoting the events is becoming more 
and more prominent with increasing amplitude range but 
narrower base. Just as in the bivariate case, the short time- 
scale fluctuations in this multivariate example can be elimi- 
nated by an appropriate (causal) moving average. Here this 
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Fig. 3. Time profile of the real-time SPIKE-distance for the two examples 
used in Fig. 2. In both subplots regular traces are depicted in blacl< while 
the causal moving averages are shown in green. A. Bivariate example. 
While the regular dissimilarity profile Sr{t) attains a local maximum 
for each spike, its subsequent decay still captures the relative spiking 
behavior. This can be seen best with the moving average which reflects 
the correct long-term trends. B. Multivariate example. The averaged 
bivariate dissimilarity profile S°{t) exhibits periods of high values in 
certain intervals, but in contrast to the original non-causal case S^(t), 
these are not spurious (see Fig. 4). 

moving average exhibits a gradual decrease, which reflects 
the consistent increase of the spike event reliability. 

It is important to note that, in contrast to the non-causal 
SPIKE-distance (see Fig. 2), the observed peaks are not 
spurious. They occur at time instants when it is not yet 
known whether there will in fact be a reliable spiking event 
or not. This is illustrated in Fig. 4 with two spiking events 
which are identical (one spike per spike train) except for 
the omitted second half of the second event. While for both 
events the instantaneous values in the first part necessarily 
have to be identical (and very high since only some of the 
neurons have recently spiked), the differences in the sec- 
ond part become evident as more and more information be- 
comes available. For the first event, once all neurons have 
fired, a rapid decrease of S'°(i) can be observed, while the 
lower firing reliability in the second event leads to a slow 
decrease. 

2.4. Practical considerations 

The improved dissimilarity profile S{t) of the SPIKE- 
distance is piecewise linear rather than piecewise constant 
as is the case for the ISTdistance. Therefore, when the lo- 
calized visualization is desired, a new value has to be cal- 
culated for each sampling point and not just once per each 
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Fig. 4. Real-time SPIKE-distance: Peaks during reliable spiking events 
are not spurious. Multivariate example with two spiking events which 
have identical first halves but different continuations. After 100 ms some 
of the neurons have just spiked while others have been silent for some 
time, and for both events this large variability is correctly reflected by 
high values of the averaged bivariate dissimilarity profile 5"(t). At these 
points begin the shaded areas which mark the spike information that is 
not yet available. Here for the first event the initial spiking in some of 
the spike trains turns into a reliable spiking event (one spike for each 
spike train) which is reflected by a rapid decrease to very low values. 
For the second event this does not happen and accordingly the decrease 
is less pronounced. 

interval in the pooled spike train. In cases where the dis- 
tance value itself is sufficient, the short computation time 
can be even further decreased by representing each inter- 
val by the value of its center and weighting it by its length. 
This is not only faster, but it actually gives the exact re- 
sult, whereas the time- resolved calculation is a very good 
approximation only for sufficiently small sampling intervals 
dt (imagine the example of a rectangular function, at some 
point any sampled representation has to cut the right an- 
gle) . The dissimilarity profile Sr (t) of the real-time SPIKE- 
distance is hyperbolic and not linear but here also the ex- 
act result can be obtained by piecewise integration over all 
intervals of the pooled spike train. 

The calculation of the SPIKE-distance consists of three 
steps: In a precalculation step, for each spike the distance to 
the nearest spike in all the other spike trains is calculated. 
Successively, for each time instant and each pair of spike 
trains, the distances of the four corner spikes are first locally 
weighted and then normalized. These latter steps involve 
matrices of the order 'number of time instants' x 'number 
of spike train pairs', which for very long datasets with many 
spike trains can lead to memory problems. The solution to 
these problems is to make the calculation sequential, i.e., 
to cut the recording interval into smaller segments, and to 
perform the averaging over all pairs of spike trains for each 
segment separately. In the end the dissimilarity profiles for 
the different segments (already averaged over pairs of spike 
trains) are concatenated, and its temporal average yields 
the distance value for the whole recording interval. 

More information on the implementation as well as the 
Matlab source code for both the ISI- and the SPIKE- 
distance (including its real-time variant) can be found at 
www.fi.isc.cnr.it/users/thomas.kreuz/sourcecode.html. 
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Fig. 5. Instantaneous clustering for artificially generated spike trains. 

A. The clustering of the spike trains changes every 500 ms (dashed 
black lines) although there are always exactly two different clusters. 

B. Matrices of pairwise instantaneous values for the four time instants 
marked by the green lines. Both the regular (top) and the real-time 
SPIKE-distance (bottom) reflect the changing cluster association in 
each 500 ms interval. Arrows indicate features described in the text. 

3. Results 



3.1. Application to artificially generated data: 
Instantaneous clustering 

By eliminating the spurious high values in S'°(t), we have 
gained reliability of for each time instant. This allows 
us to use the instantaneous matrices of pairwise spike train 
dissimilarities to divide the spike trains into clusters, i.e., 
groups of spike trains with low intra-group and high inter- 
group dissimilarity. There are no limits to the temporal 
resolution; such a clustering can, in principle, be obtained 
for each time instant. 

In Figs. 5 and 6 we show examples of instantaneous 
clustering for both the regular and the real-time SPIKE- 
distance. Both figures depict artificially generated spike 
trains which fall into different clusters. However, in contrast 
to the overall level of spike train synchrony which remains 
rather constant (results not shown), the cluster affiliation 
of the different spike trains changes every 500 ms. The spike 
trains in Fig. 5 contain four different compositions of two 
clusters, whereas in Fig. 6 the number of clusters increases 
until a state of random spiking is reached where each spike 
train forms its own cluster. As exemplified by four different 
time instants in each figure, this varying clustering struc- 
ture is correctly reflected in the pairwise distance matrices 
of both the regular and the real-time SPIKE-distance al- 
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Fig. 6. Further examples of instantaneous clustering for artificially gen- 
erated spike trains. A. Here the three spike train clusters in the first 
500 ms are followed by 500-ms-intervals of four and then eight clusters 
while the last 500 ms contain random spiking. B. Again, both the reg- 
ular (top) and the real-time SPIKE-distance (bottom) reflect changes 
in instantaneous clustering. 

though the latter only uses information from past spiking. 

Both methods can not only distinguish different clusters 
instantaneously but are also sensitive to the detailed struc- 
ture within the clusters. An example can be seen in Fig. 5 
for the 4th spike train within the second 500 ms interval 
(see arrows) . The two methods are quite similar in the first 
two columns, but they differ considerably in the third and 
fourth column. In the third column, while for the regular 
SPIKE-distance it does not matter whether the time in- 
stant is right within a spiking event or in between two events 
(compare against the first two columns), the real-time vari- 
ant clearly separates neurons depending on whether they 
have already fired or not (cf. Fig. 4). In the fourth column, 
in contrast to the regular SPIKE-distance, the real-time 
SPIKE-distance is not yet aware of the irregular cluster af- 
filiation in the last 500 ms interval. 

The main difference between the two measures is clearly 
visible in the last column where the regular SPIKE-distance 
averages over past and future behavior and thus superim- 
poses the checkered pattern of the third interval with the 
more disordered clustering of the last interval. This last in- 
terval is not yet relevant for the real-time variant, which 
only reflects the checkered pattern of the past interval. 

Similar differences can be observed in the second and 
third interval of Fig. 6. The matrices for the regular SPIKE- 
distance are quite symmetric with respect to the dissimilar- 
ities to the different clusters since the increasing distances 
between the preceding spikes are balanced by the decreas- 
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Fig. 7. Selective temporal averaging. A. Tiie spike trains are concate- 
nated from tile examples in Figs. 5 and 6. B. Eacii column (from left to 
right) depicts tiie selective temporal average over the intervals marked 
by the horizontal lines in A (from top to bottom). Arrows indicate fea- 
tures described in the text. 

ing distances to the following spikes. In contrast, the real- 
time SPIKE-distance reflects the increasing distance be- 
tween the preceding spikes only. Thus, although the spike 
time differences between adjacent groups are similar, due 
to the normalization lower dissimilarities are obtained for 
groups of spike trains whose spikes are further in the past. 

So far we have shown clustering matrices obtained at 
certain time instants. Since such a matrix exists for each 
and every time instant t, it is possible to selectively average 
over certain time intervals. These intervals do not have to 
be continuous, selective averaging over separated intervals 
is possible as well. Four examples are shown in Fig. 7, an 
average over one individual 500 ms interval, averages over 
two consecutive as well as two separated 500 ms intervals, 
and finally the average over the whole time interval. In each 
case a linear superposition of the individual matrices can be 
observed. For the whole interval (fourth column) this means 
that the four groups of 10 spike trains which frequently 
belonged to the same cluster can still be identified (see 
arrow). Due to the fact that there were two 500 ms intervals 
(the second and the fifth) where the second and the third 
spike train group formed one big cluster, these two groups 
are more closely related than the other groups, which is 
correctly reflected by the lower values of dissimilarity. 

Another option is triggered temporal averaging where 
the matrices obtained for certain trigger time instants are 
averaged. These trigger times can be obtained from external 
influences (such as the occurrence of certain features in a 
stimulus) as well as from internal conditions (such as the 



spike times of a certain spike train). An example for internal 
triggering can be seen in Fig. 8. The artificially generated 
data consist of 20 Poisson spike trains (Fig. 8A). Five of 
these spike trains were modified such that they contained 
slightly lagged and jittered copies of the spikes of the first 
spike train (which fires at a lower rate). The idea is that a 
spike in the first spike train causes a spike in these spike 
trains (with a non-constant delay) but they also have spikes 
of their own (which might have been caused by different 
spike trains) . 

These five spike trains are not easy to identify by visual 
inspection, and the overall temporal average is likewise un- 
able to do so (Fig. SB, left column). However, by using the 
spike times of the first spike train as trigger instants we are 
able to detect them. Averaging over the pairwise instanta- 
neous S"- values from all trigger instants yields the dissimi- 
larity matrix shown in Fig. SB (right column) in which an 
irregular grid of very small distance values can be recog- 
nized. 

Another representation of dissimilarity matrices are hi- 
erarchical cluster trees (dendrograms. Fig. SC). They are 
constructed as follows: First, the closest pair of spike trains 
is identified and thereby linked by a fl-shaped line, where 
the height of the connection measures the mutual distance. 
These two time series are merged into a single element, and 
the next closest pair of elements is then identified and con- 
nected. The procedure is repeated iteratively until a sin- 
gle cluster remains. The implementation of this method re- 
quires introducing the distance between a pair of clusters. 
In the single linkage algorithm used here, this distance is 
defined as the minimum over all the distances between pairs 
of spike trains in the two clusters. In the example of Fig. 
S, the dendrogram obtained from the average dissimilarity 
matrix (left) does not fall into separate clusters, whereas 
in the dendrogram of the triggered average (right) the five 
modified spike trains form one distinct cluster with the first 
spike train and can thus easily be identified. 

In the supplementary material we present a movie which 
uses the artificially generated spike trains from Figs. 5-7 
and includes instantaneous clustering, selective temporal 
averaging of individual or combined intervals, several ex- 
amples of triggered averaging as well as the corresponding 
dendrograms. As can be seen in the screenshot of the movie 
shown in Fig. 9, we added one more feature, the comparison 
of spike train groups, where the spike trains are manually 
assigned to subgroups, and a block matrix (and the cor- 
responding dendrogram) is obtained by averaging over the 
respective submatrices of the original dissimilarity matrix. 



3.2. Application to single-unit recordings from epilepsy 
patients 

In all previous examples we have used artificially gener- 
ated spike trains for which the relative levels of spike train 
synchrony were known and could serve as validation for the 
spike train distances. Here we present an exemplary appli- 
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Fig. 8. Triggered temporal averaging. A. Poisson spike trains with su- 
perimposed firing patterns: Five spike trains follow the spiking of the 
first spike train (with small amounts of jitter). Below the horizontal 
line denoting the overall average, triangles mark the spike times of the 
first spike train which are used as trigger instants. B. Dissimilarity ma- 
trices (for the regular SPIKE-distance only): Selective averaging over 
the whole interval (left) and temporal averaging over the instantaneous 
matrices at all trigger instants (right). C. Hierarchical cluster trees (den- 
drograms) obtained from the dissimilarity matrices in B. 

cation of both the SPIKE-distance and its real-time variant 
to field data for which no a-priori knowledge is available. As 
field data we chose recordings of neuronal spiking from the 
human medial temporal lobe. These recordings were per- 
formed at the University of Bonn in epilepsy patients un- 
dergoing seizure monitoring prior to epilepsy surgery. For 
a description of the data refer to Appendix B. 

In Fig. 10 we show the spike trains recorded from 42 sin- 
gle and multi-units of an epilepsy patient during an epoch 
which contained an epileptic seizure. For this particular pa- 
tient the epileptic focus was later confirmed to lie in the 
hippocampal formation of the left brain hemisphere. In this 
example the SPIKE-distance and its real-time variant ex- 
hibit a rather limited amount of variability before and after 
the epileptic seizure, while during the seizure both distances 
exhibit strong fluctuations, particularly during the second 
part of the seizure when the local field potentials exhibit 
high-amplitude rhythmic activity (Supplementary Fig. 1). 
Remarkably, both distances show a pronounced drop at the 
beginning of the seizure, at a time when only subtle changes 
are discernible in the continuous local field potential. This 
could be an indication that an increased level of synchrony 
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Fig. 9. Screenshot from a movie (see Supplementary Material) A. Ar- 
tificially generated spike trains. B. Dissimilarity matrices obtained by 
averaging over two separate time intervals for both the regular and the 
real-time SPIKE-distance as well as their averages over subgroups of 
spike trains. C. Corresponding dendrograms. 
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Fig. 10. Exemplary application of the SPIKE-distance and its real-time 
variant to single- and multi-unit recordings from an epilepsy patient 
(whose epileptic focus was in the left hippocampal formation) before, 
during, and after an epileptic seizure. The seizure onset and offset are 
marked by red vertical lines. Different brain regions (hemispheres) are 
separated by (thick) horizontal lines. LA=left amygdala; LAH=left an- 
terior hippocampus; LEC=left entorhinal cortex; LPHC=left parahip- 
pocampal cortex; and accordingly for the right (R) hemisphere. Note 
that both dissimilarity profiles start at zero due to the auxiliary spikes 
at time t = 0. 

among a population of neurons plays an important role in 
the generation of seizure activity. 

Fig. 11 shows the dissimilarity matrices for both spike 
train distances averaged over the periods before, during. 
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Fig. 11. Selective temporal averaging in field data. A. Spike train record- 
ings as shown in Fig. 10. B. Each column (from left to right) depicts 
the selective average over the time interval marked by a horizontal line 
(from top to bottom) in subplot A (before, during, and after the epileptic 
seizure). The dashed (solid) lines separate values obtained for different 
brain regions (hemispheres). 



and after the seizure. These matrices reflect the rela- 
tionships between all the neurons recorded from different 
regions of both temporal lobes. In the non-focal hemi- 
sphere, i.e. the side of the brain that does not contain 
the epileptic focus, neurons 26 to 42 exhibit a constant 
pattern of synchrony before, during, and after the seizure. 
In the focal hemisphere the general level of spike distance 
is lower during the seizure than before or afterwards, in- 
dicating an elevated level of neuronal synchrony during 
the seizure. The spike train distances between the two 
brain hemispheres (upper right and lower left quadrants) 
show a relatively high level of synchrony for the sparsely 
firing neurons (#26, 27, 28, 33, 34, 37) that is substantially 
diminished during the seizure. 

These findings are in line with standard theories on 
seizure dynamics (Engel and Pedley, 2008). The spike train 
distances thus appear to be suitable measures to describe 
the mechanisms of seizure generation, propagation, and 
termination at the neuronal level. Specifically they provide 
an opportunity to extend our mechanistic understanding 
of spatio-temporal seizure dynamics by elucidating the 
functional role of synchronization and de-synchronization 
processes. The real-time variant could furthermore be inte- 
grated into prospective algorithms aimed, e.g., at an early 
online detection of seizure occurrence (Jouny et al., 2011). 



3.3. Application to continuous data 

To our knowledge, there are no time series analyses meth- 
ods that allow to reliably track instantaneous synchrony in 
continuous data, i.e., to map their local dissimilarity to one 
value for each time instant, either once the complete seg- 
ment is available for analysis or online in real-time. Now 
that we have provided such a method for discrete data, the 
question arises whether it is possible to extend their ap- 
plicability to continuous data. Thus, in the third applica- 
tion we use the SPIKE-distance and its real-time variant to 
measure time-resolved dissimilarity in continuous signals. 
For this purpose, a high temporal resolution is once again 
beneficial since changes in synchronization can occur on 
very small time scales. Examples include the transition to 
seizure observable in the EEG of epilepsy patients (Lopes- 
DaSilvaet al, 2003). 

The principal idea is to first transform the continuous 
time series into one or more sequences of discrete events 
which are chosen to capture the most relevant characteris- 
tics of the data. In the case of neuronal membrane poten- 
tials, these are the spikes. Under the assumption that both 
the shape of the spike and the background activity carry 
minimal information, neuronal responses are reduced to a 
spike train in which the only information maintained is the 
timing of the individual spikes. For the rather oscillatory 
EEG data that we analyze here, each continuous time series 
is transformed into two separate sequences of local maxima 
and local minima; similar choices were made e.g. in Quian 
Quiroga et al. (2002) and Kugiumtzis et al. (2004). The 
SPIKE-distance is applied to both kinds of sequences, and 
the two resulting dissimilarity profiles are averaged. As be- 
fore, the temporal average of this combined profile yields 
the SPIKE-distance. For this type of data a causal calcu- 
lation is also possible, the procedure is the same as used 
above for the regular SPIKE-distance. 

To validate the instantaneous measures of synchroniza- 
tion on controlled continuous data, we used data freely 
available on the internet via the EEG time series download 
page (www.meb.uni-bonn.de/epileptologie/science/physik/ 
eegdata.html) of the Department of Epileptology at the 
University of Bonn, Germany (Andrzejak et al., 2001). We 
randomly selected = 10 independent channels (sampling 
rate 173.16 Hz) and then introduced a time-dependent 
mixing parameter m. For m — the original independent 
channels are maintained whereas ioi m ~ 1/N all channels 
are exactly identical. Intermediate values of m interpolate 
linearly between these two extremes. Although not rele- 
vant for the event detection, the variance, which due to the 
averaging decreases for small m, is adjusted to maintain 
the appearance of a regular EEG signal. 

In Fig. 12A we show an illustration of the application 
of the averaged bivariate SPIKE-distance and its real-time 
variant to continuous EEG data using just 3 channels over 
a short interval of time which includes a transition from 
independent channels to perfectly synchronous channels. 
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Fig. 12. Exemplary application of the averaged bivariate SPIKE-distance 
and its real-time variant to continuous data. A. On top we depict three 
continuous BEG signals of 0.5 s duration which exhibit a transition from 
being independent to being perfectly synchronous after about 0.5 s (the 
mixing parameter is superimposed in green). The maxima and minima 
of the BEG signals are marked in red and blue, respectively, and the 
same colors are used below for the dissimilarity profiles obtained from 
applying the two SPIKB-distances to the respective event time series 
only. For each SPIKB-distance the black trace depicts the average over 
these two profiles and the green trace its moving averages. B. On top 
we show ten mixed signals generated from 10 different BBG-channels 
by means of an oscillating mixing parameter. For maximum mixing we 
obtain identical signals, while without mixing the original BBG channels 
are preserved. Below we depict the dissimilarity profiles of the averaged 
bivariate SPIKB-distance and its real-time variant, respectively, (black) 
as well as their moving averages (green). The excerpt used in subplot 
A is marked by a small black box (see arrow). 

This transition is traced by both measures. Fig. 12B de- 
picts the dissimilarity profiles of the two SPIKE-distances 
for alternating piecewise constant variations of the mixing 
parameter which converges stepwise towards an intermedi- 
ate level. Both measures are capable of monitoring these 
transitions in the level of synchronization. Starting from 
zero values for identical synchronization and high values for 
independent channels, two gradual transitions can be ob- 
served. However, as the different gradients show, it is easier 
to trace deviations from identical signals than deviations 
from independent signals. 

4. Discussion 

There are three different contributions in this study. 
We have improved the dissimilarity profile of the SPIKE- 



distance (Kreuz et al., 2011) by eliminating the spurious 
high values that were previously obtained for event-like 
firing patterns, we have added a variant of the SPIKE- 
distance which is causal and allows real-time calculation, 
and, finally, we have extended the applicability of both of 
these variants to continuous data (such as EEG). 

The instantaneous reliability obtained by the elimina- 
tion of the spurious high values has very important con- 
sequences. In addition to the improved capability to track 
the overall level of synchronization within a group of two 
or more spike trains, it is now possible to look for in- 
stantaneous clustering, i.e., to represent evolving patterns 
of (dis)similarity in multiple spike trains either as instan- 
taneous matrices of pairwise spike train dissimilarities or 
as hierarchical cluster trees (dendrograms). The fact that 
there are no limits to the temporal resolution allows fur- 
ther analyses such as selective and triggered temporal av- 
eraging. In real multi-neuron data, triggering might help 
to detect converging or diverging patterns of firing prop- 
agation. Finally, spatial averaging over spike train groups 
is possible. In applications to real data, these groups could 
be different neuronal populations or responses to different 
stimuli, etc. 

In cases where the complete spike trains are known, 
the regular SPIKE-distance compares favorably against 
the real-time SPIKE-distance. It is locally more reliable 
because it has information about both the past and the fu- 
ture at its disposal. For the real-time SPIKE-distance the 
lack of information about spikes that occur in the future 
leads to high values during reliable spiking events. While 
these values are not spurious, they are not really informa- 
tive since they refiect local uncertainty which can easily 
be resolved by incorporating all the information available. 
However, in situations which demand an online monitor- 
ing the real-time SPIKE-distance can be relied upon as 
we could show both for simulated and field data. A good 
example for the latter is Fig. 11 in which even with only 
half the information available the results achieved were 
very similar to the ones for the regular SPIKE-distance 
(compare the upper and lower row in Fig. IIB). 

The capability to track the synchrony between two or 
more spike trains (or continuous recordings) opens up new 
possibilities for several important applications. Synchro- 
nization has been hypothesized to play a pivotal role in 
neuronal coding (Miller and Wilson, 2008), and thus a real- 
time tracker of spike train synchrony could be an essential 
tool for a rapid online decoding with brain-machine inter- 
faces in order to control prosthetic devices (Hochberg et 
al., 2006; MuUiken et al., 2008; Sanchez, 2008). In epilepsy 
patients, monitoring the spiking activity of large ensembles 
of single neurons could lead to a better understanding of 
the mechanisms of seizure generation. Furthermore, if neu- 
ronal spiking patterns turn out to be specific predictors of 
seizure occurrence as reported by a recent study (Truccolo 
et al., 2011), the real-time SPIKE-distance could be a vi- 
able tool for the implementation of a prospective seizure 
prediction algorithm (Mormann et al., 2007). 



In a similar manner the SPIKE-distance and its real-time 

variant can be applied to monitor synchrony in continu- 
ous data (which are first transformed into discrete data). 
Even the analysis of mixed continuous-discrete signals is 
possible (see also Andrzejak et al., 2011). A potential ap- 
plication could be the combined analysis of discrete spike 
trains and continuous neuronal oscillations. In particular, 
it would be interesting to investigate neuronal synchroniza- 
tion patterns in dependence of the phase of the local field 
potential, a scenario reminiscent of the 'neuronal communi- 
cation through neuronal coherence' scenario (Fries. 2005). 

Note that there is a conceptual analogy between the im- 
proved definition of the SPIKE-distance and the similarity 
measiire proposed by Hunter and Milton (2003) which basi- 
cally sums the exponentially weighted distances from each 
spike to its nearest spike in the other spike train. However, 
the main difference, apart from the additional exponential 
weighting and the parameter (the decay constant of the 
exponential) , is that the calculation is not done in a time- 
resolved manner by means of a local weighting of th(^ tc^rms 
of the four corner spikes. Instead, the term for each spike 
is considered just once. The SPIKE-distance is a rather 
elementary measure that can be regarded as complemen- 
tary to cross correlation. Both measures rely on differences 
in spike timing. However, while the latter estimates spike 
synchrony in dependence on a time lag but is not sensitive 
to instantaneous synchrony, the SPIKE-distance estimates 
instantaneous synchrony but is not sensitive to time lags. 

The SPIKE-distance and its real-time variant share 
many properties with the ISI-distance (see Appendix A) 
but there are also a few essential differences. Here we 
provide an overview of both: 

• ISI-and SPIKE-distance are complementary 

The ISI-distance is based on interspike intervals and 
quantifies covariations in the local firing rate, while the 
SPIKE-distance tracks synchrony mediated by spike tim- 
ing. This does not mean that the ISI-distance is sensi- 
tive to rate coding and the SPIKE-distance sensitive to 
temporal coding. Rather, it is the relative timing of in- 
terspike intervals and spikes, respectively, that matters. 

• Straightforward extension to multiple spike trains and 
normalization 

Both distances directly address the lack of approaches 
to analyze multiple spike train data (Brown et al., 2004). 
They can readily be applied to very large datasets (some 
of the datasets analyzed contained over 100 spike trains 
with a total of almost 250000 spikes). Both the bivariate 
and averaged bivariate distances are normalized between 
zero and one where zero is obtained for identical spike 
trains only. The same limits hold for the underlying dis- 
similarity profiles. However, there is a difference regard- 
ing the range of values obtained. For the piecewise con- 
stant ISI-distance the dissimilarity profile and the dis- 
tance itself cover the same range. In particular, the dis- 
tance value can come arbitrarily close to the maximum 



value of one. For the piecewise linear SPIKE-distance this 

is not the case. Its value is in most cases much smaller 
than the maximum instantaneous value and typically be- 
low 0.5. This is in line with the observation that for Pois- 
son spike trains of equal rate the expectation values for 
the ISI- and the SPIKE-distance are 0.5 and 0.295, re- 
spectively. 

• Different levels of information reduction 

For both distances there are three different levels of 
information reduction. The starting point is the most de- 
tailed representation in which one instantaneous value is 
obtained for each pair of spike trains. The resulting ma- 
trix of size 'number of sampled time instants' x 'squared 
number of spike trains' (i.e. jf^(tn)N^) can be appreci- 
ated best in a movie; an example can be found in the 
supplementary material. For the first step of information 
reduction there are two possible averages: the average 
over spike train pairs and the temporal average. These 
commute, either can be performed first. By averaging in- 
stantaneously over all pairs of spike trains, a dissimilarity 
profile for the whole population (e.g. Fig. 2) is obtained 
whereas the temporal averaging leads to a bivariate dis- 
tance matrix (e.g. Fig. 5). Finally, in both cases the sec- 
ond average leads to one distance value which describes 
the overall level of synchrony for a group of spike trains 
over a given time interval. Note that here we restricted 
the analysis to average values, but also the application of 
higher order statistics (such as the variance) is conceiv- 
able. 

Depending on the application in mind, the appropriate 
representation can be chosen. Mapping the similarity of a 
whole population onto one single value might allow for an 
easier comparison but one might lose too much informa- 
tion since high and low spike time differences at different 
time instants or for different pairs of spike trains might 
cancel each other out, leading to a value that could also 
be obtained for a constant intermediate level of similar- 
ity. This is one example where a higher-order statistics 
such as the variance could be useful. 

• High temporal resolution: Reliance on local information, 
then global averaging 

Both distances rely on instantaneous values which only 
take into account local information from preceding and - 
with the exception of the real-time SPIKE-distance - fol- 
lowing spikes. Temporal averaging (Eq. 1) then yields a 
more global picture by means of a single- value distance. 
The fact that the averaging window can be chosen ar- 
bitrarily allows for easy comparisons between the levels 
of spike train synchrony in different time intervals with- 
out the need of recalculations. The high temporal reso- 
lution which renders a sliding-window analysis obsolete 
compares them favorably to other widespread measures 
of spike train variability such as the Fano factor. 

• Dependence on spike train of origin 

Both distances take into account the spike train of ori- 
gin for each individual spike and neither of them is in- 



variant to shuffling spikes among the spike trains. This 

is in contrast to measures that act on the pooled spike 
train, such as all measures based on the Peri-Stimulus 
Time Histogram (PSTH), which yield the same value re- 
gardless of how spikes are distributed among the different 
spike trains and could possibly fail to distinguish between 
qualitatively different behaviors such as high reliability 
spiking and chain-fire bursting (cf. Kreuz et al., 2011). 

• Computational efflciency 

Both distances are based on simple differences and ra- 
tios which furthermore can easily be parallelized. For a 
large set of very long spike trains for which computer 
memory might be a concern, the dissimilarity profiles for 
successive intervals can be calculated sequentially. 

• Time-scale independence and absence of parameters 

Both distances are invariant to stretching and com- 
pressing of the spike trains. They are also time-scale 
adaptive since the information used at each time instant 
is not contained within a window of fixed size but within 
a time interval whose size depends on the local rates 
of the spike trains. In contrast to time-scale dependent 
measures such as the Victor-Purpura (Victor and Pur- 
pura, 1996) or the van Rossum distance (van Rossum, 
2001), parameter-free single- valued methods give a more 
objective and comparable estimate of neuronal variabil- 
ity. Other drawbacks of time-scale dependent measures 
arc the computational cost and the time and effort that 
is needed to find the right parameter. Moreover, it is not 
at all guaranteed that there exists a right parameter. 

One of the main arguments for the use of timc-scale- 
dependent measures of spike train (dis)similarity is 
their potential insight into the precision of the neuronal 
code (Victor and Purpura, 1996). This argument has 
recently been reevaluated in Chicharro et al. (2011) and 
Chicharro et al. (2012) for transient constant and time- 
varying stimuli, respectively. According to these studies 
the optimal time-scale obtained from a spike train dis- 
crimination analysis is far from being conclusive. Rather 
it results in a non-trivial way from the interplay of many 
different factors such as the distribution of the informa- 
tion contained in different parts of the response and the 
degree of redundancy between them. 

The Matlab source codes for the ISI- and the SPIKE- 
distance (including its real-time variant) ) as well as 
additional material (including the supplementary fig- 
ure and the supplementary movie) can be found at 
www.fi.isc.cnr.it/users/thomas.kreuz/sourcecode .html. 

Appendix A. The ISI-distance 

While the dissimilarity profile of the SPIKE-distance is 
extracted from differences between the spike times of the 
two spike trains, the dissimilarity profile of the ISI-distance 
(Kreuz et al., 2007, 2009) is calculated as the instantaneous 
ratio between the interspike intervals ccjgl and a;|g| (Eq. 5) 



according to: 

= ^lsl(^)/^isl(J) - 1 if < 4Uit) 

\ — {xi^i{t)/xjgj{t) — 1) otherwise. 

(A.l) 

This ISI-ratio equals for identical ISI in the two spike 
trains, and approaches —1 and 1, respectively, if the first 
or the second spike train is much faster than the other. 
For the ISI-distance the temporal averaging analogous to 
Eq. 1 is performed on the absolute value of the ISI-ratio, 
thus both kinds of deviations are treated equally. Since the 
ISI-distance relies on the instantaneous ISI- values and thus 
requires knowledge about the following spikes, no causal 
real-time extension is possible. The dissimilarity profile is 
condensed into a distance value by means of temporal av- 
eraging analogous to Eq. 1. In Andrzejak et al. (2011) the 
ISI-distance has been integrated in a measure which can de- 
tect unidirectional coupling not only between spike trains 
but also between spike trains and time-continuous flows. 

Appendix B. Field data: Single- unit recordings 
from epilepsy patients 

All experimental recordings were performed prior to and 
independently from the design of this study and conformed 
to the guidelines of the Medical Institutional Review Board 
at the University of Bonn. Electrode locations were based 
exclusively on clinical criteria and were verified by MRI 
and by computer tomography co-registered to preoperative 
MRI. Each electrode probe had eight high-impedance (typ- 
ically 800 to 1000 kil) micro-wires (Platinum/Iridium, 40 
/im diameter) protruding from its tip (Adtech Inc, Racine, 
MN) . The differential signal from bipolar montages of the 
micro- wires (one wire was used as local ground) was ampli- 
fied using a 128-channel Neuralynx system (Digital Lynx 
lOS, Neuralynx, Bozeman, MT), filtered between 0.1 and 
9000 Hz, and sampled continuously at 32 kHz for later pro- 
cessing and analysis. The Neuralynx head stages used had 
unity gain, very high input impedances (ca. 1 Tfl) and no 
phase shift. Spike detection and sorting was performed af- 
ter band-pass filtering the signals between 300 and 3000 Hz 
(Quian Quiroga et al., 2004). 
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Caption Supplementary Figure: An epileptic seizure 

recorded with microwires located in different regions of the 
medial temporal lobe. Some recording channels (e.g., LAS) 
are missing because wires were broken over the course 
of the EEG-vidco monitoring. Red vertical lines denote 
seizure onset. Seizure activity remained confined to the 
focal (left) hemisphere. Note that the neuronal spiking 
activity shown in Figs. 10 and 11 was extracted from the 
same signals after highpass filtering. LA=left amygdala; 
LAH=left anterior hippocampus; LEC=left entorhinal cor- 
tex; LPHC=lcft parahippocampal cortex; and accordingly 
for the right (R) hemisphere. 



Caption Supplementary Movie: This movie uses the ar- 
tificially generated spike trains form Figs. 5 - 7. In the first 
part the instantaneous dissimilarities and the clustering 
dendrograms of both measures are updated as the green 
time bar moves from left to right. The second part com- 
prises several examples of selective temporal averaging of 
individual or combined intervals as well as triggered aver- 
aging. A. Artificially generated spike trains. B. Dissimilar- 
ity matrices obtained by averaging over two separate time 
intervals for both the regular and the real-time SPIKE- 
distance as well as their averages over subgroups of spike 
trains. C. Corresponding dendrograms. 
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